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, Abstract 

(N . 

It is shown that natural boundary conditions for non-relativistic wave 
functions are of periodic or of homogeneous Robin type. Using asymp- 
totic central symmetry of Hamiltonian and theory of singular differential 
equations the many-electron wave function is expanded in series both 
in the vicinity of Coulomb singularities and at infinity. Hydrogenic an- 
gular dependence of three leading terms of expansion about Coulomb 
. singularities is found. Exact first- and second-order cusp conditions are 

O |' obtained demonstrating redundancy of spherical average in Kato's cusp 

condition. Our first-order cusp condition exhibits CP symmetry. Ho- 
C"j . mogeneous Robin boundary conditions are obtained for aperiodic many- 

electron systems from the expansions. Use of our explicit boundary con- 
ditions improves both speed and accuracy of numerical calculations. A 
confluent hypergeometric series defining arbitrarily high order cusp con- 
ditions for the spherically averaged Hamiltonian is presented. 

> 

O 1 Introduction 
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O 

Boundary conditions play important role in eigenvalue problems of mathemati- 
cal physics even if they are imposed implicitly. Explicit use of boundary condi- 
tions is crucial to numerical calculations. The role of regularity and boundary 
conditions in existence of quantum mechanical eigenvalue problems was first 
recognized by Schrodinger [1] and von Neumann [2]. Boundary and regularity 
conditions for the Schrodinger equation were first studied in detail by Jaffe [3J H] 
who recognized that both ip and its first derivatives should be bounded even at 
singular points of the potential. McCrea and Newing [5] discussed boundary 
5^ ■ conditions for the if^" molecular ion. 

In numerical calculations of many-electron systems, explicit boundary con- 
ditions should be imposed on the wave function at each singular point of the 
Hamiltonian, namely, at Coulomb singularities and at infinity The exact many- 
electron wave function satisfies boundary conditions not only at the nuclei but 
also at electron-electron coalescence points r.^ =0 (i ^ j). 

Asymptotic behaviour of many-electron wave functions was discussed in the 
framework of Hartree-Fock (HF) approximation by Handy, Marron and Silver- 
stone [5] and in the general case by Katriel and Davidson [7] and by Patil [HI [H] ■ 
Vitanov and Panev [TU] extended the investigation to excited atomic states. The 
long-range behaviour of density was studied by Tal [TT] and by Levy, Perdew 
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and Sahni [15] . Ernzerhof, Burke and Perdew [13] discussed the long-range be- 
haviour of ground-state wave functions, one-electron density matrices, and pair 
densities including their angular dependence. 

The study of many-electron wave functions about Coulomb singularities 
dates back to Lowdin [14] who first recognized their cusp-like behaviour. He 
obtained, by analysing tabulated atomic HF data, hydrogenic r-dependence 
Rni(r) = r l j n i(r) oc r e [l — Zr/(l + 1) + . . .] of radial wave functions in the 
vicinity of nucleus resulting in cusp condition 



fn j 
fnl 



'+1 



(1) 



Similar conditions were obtained both for electron-nucleus and electron-electron 
coalescences in the ground state of the He atom by Roothaan and Weiss [TS] by 
analysing singularities of Hamiltonian which were used as constraints on trial 
functions in their variational calculation in order to avoid divergence of local 
energy Ei oc = ^~ 1 H i ^ at coalescence points r\ = 0, ri = and ri2 = 0. They 
have introduced the reduced mass of coalescent particles into the cusp condition 
(see quotation from p. 196 of Ref. [T5]): 

Each cusp refers to a pair of Coulombic particles; the sign of the 
cusp is positive for a repulsion, negative for an attraction, and the 
magnitude is equal to the product of the two charges and the reduced 
mass of the two particles. 

Roothaan, Sachs and Weiss [16] used above cusp condition in variational calcu- 
lation of light atoms and ions (up to Z = 10). The first variational molecular 
calculation using the cusp constraint (for H2) was performed by Kolos and 
Roothaan [17] , 

The near-nucleus behaviour of an arbitrary many-electron eigenfunction ^ = 
^(r|r2, . . . , rjv) was given by Kato's [TS] theorem lib 



Or 



■Z9(r = 0), (2) 



where Z denotes nuclear charge of the atom located at the origin and overline 
symbol stands for the spherical average. Steiner [Tj5] added a corollary to the 
theorem which relates the charge density and its derivative at the nucleus of an 
atom in a similar manner. M. Hoffmann-Ostenhof and Seiler I20J generalized 
Kato's result to multiple coalescences; their fixed-nucleus approximation was 
removed by Johnson [51] using Jacobi coordinates. Validity of Eq. ([5J is limited 
to the s states since wave function ip oc r l has a root of multiplicity I at the 
origin. Unfortunately, Eq. ((5J and its generalizations cannot be used as exact 
boundary conditions due to the spherical average. It will be shown in section [4] 
that the spherical average is redundant in Kato 's cusp condition ([2]) which was 
already indicated for above mentioned few-electron systems by Roothaan and 
his co-workers [151 HU Q7] . 

Results of our thorough and rigorous analysis differ from that of three highly- 
cited papers [22] [23] [24] but are in accordance with that of Lowdin [M] and of 
Roothaan and his co-workers [TSJ [TH] HZ] . Pack and Brown [55] added the phrase 
"all /,m" to Lowdin's cusp condition ([1} by solving many-electron Schrodinger 
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equation in the vicinity of coalescence points which is not justified by our more 
rigorous solution. Bingel |23l 124] introduced an additional linear term propor- 
tional to the cosine of an unknown angle in the series expansion of wave function 
of which spherical average is always zero. In our rigorously derived expansion, all 
three leading terms exhibit the same hydrogenic angular dependence Yi m (i), ip) 
as a result of isotropy of three leading terms of many-electron Hamiltonian in 
the vicinity of Coulomb singularities. Perhaps some recent research (e.g. Refs. 
[251 126] ) relying on these old results should be revised. 

As a result of spherical average introduced by Kato the higher-order cusp 
conditions were derived only for the spherically averaged Hamiltonian by Ras- 
solov and Chipman [37] which were generalized to the excited states by Nagy 
and Sen [28, 29, 30 . The near- nucleus anisotropy of Hamiltonian was first taken 
into consideration by Qian and Sahni [3T], nevertheless omitting term O (r°). 
We will show in sections 13.31 and 14.21 that first and second order cusp condi- 
tions are exact without spherical averaging and only that of third and of higher 
orders need spherical average. In addition, we present a confluent hypergeomet- 
ric series in section 13.41 defining arbitrarily high order cusp conditions for the 
spherically averaged Hamiltonian. 

Fournais, M. and T. Hoffmann-Ostenhof and S0rensen [32] recently decom- 
posed the many-electron wave function as V'( x ) — ^^'(x) + |^i| ip( 2 \x) in the 
neighbourhood of Coulomb singularities using inconvenient Cartesian coordi- 
nates x = (xi, . . . ,xn) € M. 3N inhibiting the recognition of asymptotic central 
symmetry of Hamiltonian. In addition, their component functions an d 
are real-valued. Our decomposition (|2"3"|) is slightly different due to the isotropy 
of singularity taken into consideration: the first term is an eigenfunction of the 
spherically averaged Hamiltonian, the second term reflects the anisotropy of 
molecular or crystalline potential. 

Since Kato's rigorously proved spherically averaged cusp condition @ con- 
tradicts Eq. ([T]) which was successfully used in few-electron calculations by 
Roothaan and his co-workers, a different but similarly rigorous treatment is 
needed to decide between the two forms. The primary aim of this paper is to 
derive boundary conditions by solving many-electron Schrodinger equation both 
at large distances and in the neighbourhood of Coulomb singularities using a 
well-tried formalism of mathematical physics, the theory of singular differential 
equations. The key idea of this paper is that due to asymptotic isotropy Qof 
many-electron Hamiltonian both at Coulomb singularities and as r — > oo the 
angular momentum asymptotically commutes with the Hamiltonian leading to 
asymptotically hydrogen-like wave functions with definite values of quantum 
numbers I and m. 

The outline of the paper is as follows. In section ® the allowed boundary 
conditions of non-relativistic one-particle quantum mechanics are studied. In 
section 13 - 1 1 the isotropy of singularities of many-electron Hamiltonian is dis- 
cussed. In section 13.21 the many-electron wave function is expanded at large 
distances in terms of irregular solid spherical harmonics (]13j) . At the beginning 
of section 13.31 the many-electron Hamiltonian (Ti"4"]) is studied in the vicinity of 
Coulomb singularities, where explicit form (I15j) of potential contribution of non- 

1 Term asymptotic equality is used in this paper in the sense of liniz-^Q [f(%)/g{%)] = 1. 
Isotropic and anisotropic terms of expansions are distinguished throughout this paper so 
a boldface argument of ordo symbol O (r n ) indicates the anisotropy of omitted terms of 
expansions. 
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coalescent particles is obtained. Then it is shown that due to isotropy of leading 
terms of Hamiltonian, all three leading terms of many-electron wave function 
exhibit hydrogenic angular dependence (|23[) . therefore first and second order 
cusp conditions (|21al I21b|) are exact without spherical averaging. The section 
ends with the decomposition of many-electron wave function as the sum of two 
functions, where the first term is an eigenfunction of the spherically averaged 
Hamiltonian, the second term reflects the anisotropy of molecular or crystalline 
potential. In section [3T4l eigenfunctions of the spherically averaged Hamiltonian 
are discussed leading to recurrence relations (|26al I26bl) defining cusp conditions 
of arbitrarily high orders for this special case. In section 14.21 exact first and 
second order cusp conditions (|30l I31al I31bl l32j) are presented demonstrating 
the redundancy of spherical average in Kato's cusp condition @. As a result 
of asymptotic hydrogenic behaviour of many-electron wave function both at 
Coulomb singularities and as r — > oo we obtain boundary conditions of homoge- 
neous Robin type (|2"51 |2T))) similarly to Eqs. © and §§§ derived in section [5] for 
the H atom. In section [5j physical and numerical consequences of our results 
are discussed. 

2 Boundary conditions in one-particle quantum 
mechanics 

The problem of hydrogen atom was solved by Schrodinger in his historically 
famous paper pQ without imposing explicit boundary conditions on the wave 
function. In the first version of manuscript, the requirement of stationarity 
of current flux was used as a constraint on the variational problem which was 
changed to the weaker normalization condition by an addendum (cf. equations 6 
and 24 of Ref. [1]). Latter form is more conventional mathematically since it is 
compatible with Sturm-Liouville theory of eigenvalue equations, where J \ip\ dv 
represents the denominator of Rayleigh quotient. Von Neumann has concluded 
that normalization condition for the wave function and requirement of self- 
adjointncss of the Hamiltonian is equivalent to imposing both boundary and 
regularity conditions on the wave function [2J. 

In fact, these conditions are too weak to enforce unique regular solutions of 
Schrodinger equation. The normalization condition does not exclude irregular 
particular solution ijj cx r~ l ~ x for s states of the Coulomb problem [33J hence it 
is excluded by hand both in Schrodinger's paper [T) and in the textbooks. The 
requirement of self-adjointness does not lead to a unique eigenvalue problem 
since a little-known theorem (341 135] of Sturm-Liouville theory of differential 
equations states that any of following two types of boundary conditions are 
consistent with self-adjointness of the Liouville operator: 

1. periodic boundary conditions 



^(a)-i>(b) = 0, 



(3a) 
(3b) 



2. homogeneous Robin boundary conditions 

aiij)'(a) + /3iip(a) = 0, 
a 2 ^{b)+M{b) = 0, 



(4a) 
(4b) 
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where a's and /3's are real constants, a and b denote endpoints of the interval 
0. The theorem can be generalized to partial Sturm-Liouville equations by 
taking function values and normal derivatives over hypersurfaces of the domain 
(the proof is based on the definition of self-adjointness and the use of Green's 
theorem). 

Equations (f5a1l3b|) are known as Born - von Karman [SHUSHIHE] or Bloch [33] 
boundary conditions of solid state physics. In view of above and Bloch's theo- 
rem we can state that eigenf unctions of aperiodic systems satisfy homogeneous 
Robin boundary conditions. Homogeneous Dirichlet and Neumann boundary 
conditions for model problems of the textbooks are special cases of (Hal I4b|) . 
Boundary conditions (Hal I4bl) can be divided by arbitrary constants so coeffi- 
cients a j / y/af + 01 = sin 7$ and ft / y/af + (3f = cos 7, define angles 71 and 72 
representing the boundaries. Random coefficients for amorphous materials, in 
the Wannier [40] representation, may be interpreted as random walk of these 
"phase points" around the unit circle which leads to a band structure similarly 
to the periodic boundary conditions. 

As an example of Eqs. (|4al I4b[) let us recover hidden boundary conditions 
for a non-relativistic £f-like ion with nuclear charge Z using known properties 
of hydrogenic bound-state wave functions. The normalization condition guar- 
antees a part of boundary conditions, namely, vanishing at infinity. Asymptotic 
exponential decay of wave function tp — tp(r) is described by limit of logarithmic 
derivative 

I dip 



lim = -V=2E, (5) 

r->oo 1p Or 

where, and throughout this paper, atomic units (h = e = m e = Aire® = 1) are 
used. The normalization condition results in a mild singularity (hyperconical 
cusp) of the wave function: whereas it is continuous everywhere, its directional 
derivatives are bounded but discontinuous at the Coulomb singularity. Since 
eigenfunctions of the central field problem are separable as ip — Rnii^Yem^, tp), 
their directional logarithmic derivatives are of the form 

e • Vi/j e • e r dR n £ e ■ e$ dYi m e • e v dYi m 



ip R n t dr Y em r dd Y tm r sin dip 

where e stands for unit vector of the selected direction, {e r , e^,e v } is the basis of 
the spherical polar coordinate system. In radial directions defined by e-e r = ±1 
and e ■ e§ = e ■ e v = 0, above expression reduces to 

e ■ W ±1 dR ne ±1 dtp 



tp R n i dr ip dr 

Since radial wave function R n g = r l u n iir) of the central field problem has a 
root of multiplicity I at the origin the l'Hospital rule should be applied t times 



2 The natural boundary condition for the momentum operator — i^jr, discussed by von 
Neumman [2] as an example, is simply J3a.l l . Widely used substitution ip = u(r)/r eliminates 
not only the first derivative from the Laplacian but also the periodic boundary conditions 
since r equals the square root of coefficient p = r' 2 of Liouville operator ^ [p( r )^j "t" 9( r )> 
i.e. only one branch of i^/p is used which illustrates the importance of single-valuedness of 
coefficients in the theory of differential equations. In d-dimensional hyperspherical coordinates, 
the substitution eliminating the first derivative is of the form tp = u{r)/A d ~ 1 ^ 2 , where the 
square root appears explicitly. 
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in order to obtain a definite limit 

ltaJ-^ = lJm^ = (/ + l)^ = -3 
r^oR ni dr r^O R W V 'u ne (0) 

where differentiation rules (I35al I35b|) and explicit form of the hydrogenic radial 
wave function [f| are used. The discontinuity of the radial logarithmic derivative 
at the nucleus is then characterized by 

1 dtp d e+1 ip 
lim — — = lim r = (6) 
r^±0 ip dr r^±0 df.ip 

where 

r -> +0 = (r -> 0, t?, y>), 

r -> -0 = (r -> 0, 7r - 7r + (ys). 

We emphasize that all steps of above derivation, except the last one, rely on 
isotropy of singularity of Hamiltonian at r = 0. It is interesting to observe 
the CP invariance of cusp relation ©. Pair of ([5} and © obviously represent 
homogeneous Robin boundary conditions of the form (|4al I4b[) . 

Similar boundary conditions will be obtained for many-electron wave func- 
tions in section [4] as a result of asymptotic central symmetry of many-electron 
Hamiltonian both in the vicinity of nuclei and at large distances. 

3 Behaviour of many-electron wave function at 
singular points of Hamiltonian 

Let us consider non-relativistic Hamiltonian describing N particles interacting 
with each other by Coulomb potentials 

N * N-l N 
i— 1 i— 1 j— z+1 J 

where rrij and g.; denote mass and charge of the z-th particle, respectively (m^ = 1 
and qi = —I for electrons and m t = 1836^^, % = Zj,, i/ = 1, 2, . . . < N for 
nuclei). Electrons and nuclei will be distinguished only in the final results. 
Spin coordinates are omitted for simplicity. Many-particle wave functions \& = 
(ri, T2, . . . , rw) satisfy stationary-state Schrodinger equation 

= (8) 

which is singular both at coalescence points r = and as r — > oo. Singu- 
larities of the Coulomb potential are removable and singularity of the kinetic 
energy at infinity is essential. The singularities and other leading terms of the 
Hamiltonian are isotropic in both limiting cases. In the vicinity of coalescence 
points, the singular potential contribution of coalescent particles O (r -1 ) and 
bounded leading term of potential due to remaining non-coalescent particles 
O (r°) are isotropic. At infinity, only Coulombic monopole term O (r _1 ) is 
isotropic, whereas multipole terms O (r -2 ) are anisotropic. 



3 The power series expansion of the hydrogenic radial wave function about the nucleus is 
1 



of the form R n i = C n ( r l 



Zr i 2n^+l + l / Zr_ 
t+1 " 1 " 2(£+l)(2£+3) V ™ 
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3.1 Isotropy of singularities 



Particular solutions of a linear ordinary differential equation in the neighbour- 
hood of an isolated singular point Xq = are of product form y — f(x)g(x), 
where f(x) ensures the correct behaviour of y at the singularity and g(x) is a 
single- valued analytic function (or has at most a logarithmic singularity) which 
is non-zero at the singular point. Function f(x) is typically power, exponential 
or Gaussian function. For isotropic singularities of a linear partial differential 
equation, the solution has the form y = f(r)g(r) with r € K n and r — |r|, where 
/(r) is responsible for correct behaviour of y about the singular point and g(r) 
(reflecting anisotropy of the coefficient functions) is non-zero at the singular 
point. (One may expand it in 3 dimensions in terms of regular r x Y\^(-d, ip) or 
irregular r~ x ~ 1 Y\ )Ji (i!), ip) solid spherical harmonics.) It will be shown in this sec- 
tion that due to isotropy of singularities of Hamiltonian ([7]) , the wave function 
about singular points has the following limiting forms 

* — > r^(r), Ui (0)^0, i=l,2,...,N, 

r— >ri— 

* — > e ar r p v(r), lim v(r) ^ 0, a < 0. 

r— >oo r— >oo 

In many-body systems, in contrast to the central-field problem, only energy 
E is conserved throughout the configuration space, square L 2 and projection L z 
of angular momentum are conserved only at the singularities 

[H, L 2 ] — > 0, [H, L z ] — > 0, (9a) 

[H, L 2 ] — >■ 0, [H, L z ] — >. 0, (9b) 

where Hamiltonian ([7]) is rotationally invariant. In the neighbourhood of sin- 
gular points angular momentum quantum number I and magnetic quantum 
number m have definite values. Eigcnfunctions of many-electron Hamiltonian 
approach eigenfunctions of angular momentum when approaching singularities 
of the Hamiltonian leading to asymptotic hydrogenic angular dependence 

* — ► RitfYmntfatp), (10a) 

r— >ri— 

* — > i?ocWr o (fob) 

r— >oc 

of the wave function, where the power and exponential functions reflecting the 
singularities are included in the radial functions. In other words, molecular 
symmetries manifest only at molecular distances, where many-body Hamilto- 
nian does not commute with angular momentum. 

Due to asymptotic TJ-like behaviour of many-electron wave functions the 
equations of this section and their solution methods are similar to the mathe- 
matically correct original treatment of H atom by Schrodinger and Weyl (see 
acknowledgement in footnote on p. 363 of Ref. [I]) based on theory of singu- 
lar differential equations which differs from the simplistic textbook-derivations. 
The Coulomb singularities are treated using Fuchs' theorem exactly the same 
way as in Ref. pQ. The singularity at oo is treated using Hamburger's theorem 
which is shorter than Schrodinger's solution based on the Laplace transform. 
Of course, our equations lead to wave functions of the H atom in the special 
case of two Coulombic particles. 
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3.2 Asymptotic behaviour at infinity 

Let us consider an electron, say particle 1, separated from the rest of the system 
(mi = 1, q\ = —1 and r\ > r 2 , ■ ■ ■ , rjv). Let us introduce reduced mass M 1 and 
center of mass R' of the whole system, center of mass R" of particles except 
the electron located at ri and their separation r: 



S77 = Z^— ' R = " ' R 



r = r i R" = (r,i3,ip) = (r, w) . 
Using Laplace expansion for r > r$ 



1 OO A . ^ 

E E otzt^t^c^m 



|r-r s | ^ ^ 2A + lr A + 
the potential energy of the system can be expressed at large distances as 

TTl \ ST 47r *! ^Hv^ Xv i \ 

U ^ = E 2 A + 1 r^+i Lft r i^("i) 

A,/i i=2 

- ^I> + °(r- 2 H-^ + 0(r-), 

i=2 

where Q = J3j=i 9« by noting that Q = for neutral systems. 

Due to asymptotic isotropy f|9b[) of Hamiltonian the many-particle Schrodinger 
equation is asymptotically separable (|10bp in terms of spherical polar coordi- 
nates as r — > oo. Asymptotic radial wave function i? = i? (r|r2, ■ ■ • , rjv) satisfies 
differential equation 3 



2M' 



«±±-JS + Ofr- 



R — ► 0, (11) 



which is not an eigenvalue equation and E is eigenvalue of ©. The equation 
has an isolated essential singularity at infinity since transformation of variable 
z = 1/r leads to second order differential equation 

z 4 R"{z) + 2M' [E+(Q + l)z + O (z 2 )] R(z) —>0 

with an isolated fourth order pole at the origin. The pole being independent 
of potential is a consequence of the Laplacian. Radial equation (|TTj) is of Ham- 
burger type [4T1 1421135] 

d 2 R ( a± ao \ dR / b\ b 2 

■ a + — + 4 + ...)— +[6o + J: + 4 + ...LB=0 



dr 2 V r r 2 / dr V r r 2 
with 

at = 2, a Q =a 2 = a 3 = ... = 0, b = 2M'E, b 1 =2M'(Q + l), 



4 Centrifugal kinetic energy term is omitted since it has the same order as omitted 
anisotropic dipole term of the potential energy. 
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&2 and higher order coefficients are anisotropic. Due to essential singularity of 
the equation the solution vanishes transcendentally as r — > oo. Seeking the 
solution in the form 

R = e ar u(r), lim u(r) ^ 



we obtain 

d 2 u 
dr 2 



1\ du 
2 I a+- — - 

r / dr 



2 , 2a + bi . _ 2N 
a 2 + 6 + — + O (r 2 ) 



it — > 0. 

r— ^oo 



By equating leading term of coefficient of u{r) to zero we obtain indicial equation 
of which roots are 

a = ± v / -&o, 

where only the minus sign leads to a bounded solution. Non-essential singularity 
of the latter differential equation can be removed by substitution 

u = r^v(r), lim v(r) =/= 



yielding 



d 2 v / /3 + l\ dv 
+ 2 a + — 1 



dr 2 



dr 



2a(/3 +l) + &i 



+ O (r- 2 ) 



w — >■ 0. 

r— >oo 



By equating leading term of coefficient of v(r) to zero we obtain indicial equation 
with only root 

2a 2y/^ 

In case of a central-symmetric problem the solution of latter differential equation 
would be of the form 



v = v + vir 1 + v 2 r 2 + 



(«o 0) 



but since hi is anisotropic we restrict our solution to vq. The radial wave function 
at large distances is then 



M'(Q + 1) 



R 



(12) 



This behaviour is in accordance with results of Katriel and Davidson [7] who 
derived it less rigorously in two different ways. Since above function completely 
characterizes the singularity of Hamiltonian at oo the many-electron wave func- 
tion can be expanded in terms of irregular solid spherical harmonics as 



if = c -V-2M'Er r ^ 



! ' E ^ 



Y Xll (0,ip), (13) 



where i>oo = \/4ttvo — constant and coefficients (r~ 1 \v2 1 . . . , r^r) are single- 
valued analytic functions of r . Molecular or crystalline symmetries are re- 
flected by relations between coefficient functions (r^ 1 ^, . . . r^) for A > 0. 
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3.3 Local behaviour at Coulomb singularities 

Let us focus now our attention on coalescence of any two particles, say 1 and 2, 
while keeping remaining particles separated from them: n, r 2 < r^, . . . ,rjy 
In order to explore symmetry properties of the Hamiltonian it is convenient to 
use Jacobi coordinates by introducing reduced mass M, center of mass R and 
separation r of these two particles: 

M _ 171x1712 R = mi£l + 1712*2 

mi + m 2 ' mi + to 2 ' 

r = n -r 2 = (r,$,(p) = (r,cj). 

Many-body Hamiltonian ([7]) can be partitioned as 



where 



i=3 v 



A N A W 

2 (m 1 + m 2 ) ^2™, ^ r tf ' 



Use of Laplace expansion for r < r' 



yields 



I Ajr r"^ 

t^ = E E ^tt^ti^H^^) 



^ = E £tt [<^;m + - 

N 

= (r 3 , ...,rjy) = ^ -£p-Y^(u>; 



Expressing ri and r 2 with Jacobi coordinates r and R, putting origin of the 
coordinate system to center of mass R, using inversion property 

Y Xfi (ir-#, 7r + p) = (-l) x Y x ^ 1 cp) 

and addition theorem 

Afjr 

£ n M (^,^)n;(^^) = ^(cos 7 ,), 

/i— — A 

cos 7i = cos cos i?j + sin i3 sin i?, cos((^ — ^ ) 



5 Other authors e.g. of Refs. |22l 1271 assume well separated particles ri, r-2 <C r3, . . . ,rjv 
which is unnecessary. 
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we obtain 



W 



E 

A=0 



rn 



(Mi 



N 

\A \ "* 



2=3 



T P A (cos7i 



W (r 3 , . . . , tat) + Wi(r/r, r 3 , . . . , i-jv) r + . . . 



(15) 



for potential energy contribution of non-coalescent particles. For identical co- 
alescent particles, all odd powers of separation vanish in the vicinity of the 
coalescence point due to the inversion symmetry. This property has a profound 
consequence in the behaviour of electron-electron potentials which will be dis- 
cussed in a separate paper. Leading term 



N 

W = W (r 3 , ...,r N ) = (qi+ q 2 ) ^ 

i=3 



(16) 



of expansion (TT5|) depends on distances rather than positions of the non-coalescent 
particles. In practical calculations, the terms of above sum should be evaluated 
as expectation values |<7i/rj| \&), however the fixed- nucleus approximation 
leaves them in their original form for the nuclei. Since average of Legendre 
polynomials with A > vanish the spherical average of W equals the leading 
term of its expansion 

W = Wo (17) 

within the convergence radius of expansion (TT5"j) . Consequences of this property 
will be discussed in detail in section 13.41 Figure [T] illustrates the behaviour of 
potential W in the vicinity of a Coulomb singularity. 

Since term W of Hamiltonian (Ti"4"|) acts on separation r of particles 1 and 
2 and term G acts on their center of mass R the many-particle Schrodinger 
equation is separable resulting in an effective one-body problem. The wave 
equation for ^ = ^ (r|r3, . . . , rjv) approaches homogeneous differential equation 



A 

2M 



9i 92 



-^77 + — + W -E + O(r) 



•0 



in the vicinity of coalescence points, where E is an eigenvalue of many-electron 
Schrodinger equation ©. Due to local isotropy (l^a|) of Hamiltonian (TT4"]) the 
equation is locally separable in terms of spherical polar coordinates and the 
wave function exhibits local hydrogenic angular dependence of the form (|10a[) 
about the origin. The corresponding radial wave equation 



-1 
2Mr 2 



<7i<72 



Wo - E + O (r) \ R 



•0 



has an isolated singular point at the origin. Frobenius normal form of this 
equation is 



- 2 R" + rP{r)R' + Q(r)R = 0, 



where 



P = 2, 



-1(1 + 1) - 2M qi q 2 r - 2M(W Q - E)r 2 + O (r 3 ) 
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W =W(x=0) 














1 , 1 





-2-1012 



Figure 1: Illustrative example of a linear chain of homo- nuclear atoms showing 
the behaviour of potential term W (fT5j) of Hamiltonian (ITi|) in the vicinity of 
an electron-nucleus coalescence point x = 0. 



The singular point is removable since -P(r) and are single- valued analytic 
functions. Fuchs' theorem [33J [351 S3] states that in the neighbourhood of re- 
movable singularities, the fundamental system of solutions is 

Ri = r Xl u(r), 

i?2 = r A2 [v(r) + a u(r) lnr] , 

where SRAi > 3?A 2 , u(r'), v(r) are single-valued analytic functions, u(0) ^ 0, 
f (0) 7^ and a is a constant @. The theorem distinguishes three cases for 
existence of logarithmic term of the second solution depending on difference 
Ai — A 2 . In order to determine exponents Ai and A2 we seek the solution in the 
form r x u(r). The substitution gives indicial equation 

[A(A - 1) + 2A - iU + 1)] u(r) + O (r) — ► 

1 — >o 

of which roots are Ai = i and A2 = — I — 1. Since difference Ai — A2 = 2£ + 1 
is a non-zero integer there is no general rule for existence of logarithmic term 
hence value of a should be determined individually by substituting R2 into the 
differential equation leading to 

2M qi q 2 v(r) + 2lv'(r) + O (r) 

ol — — y — v ■ — y U 

r^o (2£ + l)u(r) + O (r) r^o ' 

i.e. the logarithmic term of the second solution vanishes at the origin hence the 
fundamental solutions are simply R\ = r e u{r) and R2 — r~ l ~ 1 v(r). Since R2 is 

6 Transformation of r to a dimensionless variable required by the logarithmic term is omitted 
for brevity. 
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unbounded at r = the physical solution is 



R—>r*u(r). (18) 

Substituting this expression into the radial equation we obtain differential equa- 
tion 



21 + 2 . 
— —v! - 2M 



qiQ2 ' W -E + O(r) 



u — >0 (19) 



r-s-0 



r 

of which solution should be analytic according to Fuchs' theorem hence it can 
be expanded in power series 

u(r) = u(0) + u'(0)r+^^r 2 + u(0) ^ 0. (20) 

Inserting it into the differential equation we obtain algebraic equation 

-2 [Mq iq2 u(0) -{£ + 1)«'(0)] r" 1 
+ [2M(£ - VKo)m(O) - 2M qi q 2 u'(0) + (2£ + 3)tt"(0)l r° + O (r) — >■ 

which can be satisfied only if 

_ u'(0) Mq iq2 

a = — — — = , 21a 

u(0) £+1 ' 1 ' 

W(0) _ (£ + l)a 2 + M(W -E) 

b ~ 2 u(0) " 2£ + 3 ' [ZLb) 

where cusp condition (|21al) removes singularity O of the equation and 

relation (|21b[) represents a constraint on curvature of the wave function at the 
origin 0. We have to note that b = b (r 3 , . . . , rjy) depends only on charges and 
distances of the non-coalescent particles similarly to Wq = Wo (rs, . . . , rjf). 

Since singularity of many-particle Hamiltonian (I14j) is isotropic at r = its 
eigenfunctions have the form 



* = r e u (r|r 3 , ...,r N ) = r e J2 r X u x „(r)Y x ^, <p) (22) 

A.// 



in the vicinity of coalescence points which can be considered as the spatial 
generalization of Frobenius series |f| , where u\^(r) = u\^(r\r3, . . . , rjy). In view 
of local commutativity (|9ap of angular momentum with Hamiltonian both £ 
and m have definite values at the origin. Since many-particle wave function is 
antisymmetric under interchange of any two electrons, £ takes only even values 
for relative singlet spin states si + s 2 = and odd values for relative triplet 
spin states si + s 2 = 1 of the coalescent electrons. Since leading terms £{£ + 
l)/2Mr 2 , qiq 2 /r and Wq of Hamiltonian are isotropic, three leading terms of 
above expansion exhibit hydrogenic angular dependence 

* — > r e u e (0) (1 + ar + br 2 ) Y lm {$, <p) + . . . . 

? >0 



^Leading expansion coefficients of radial wave functions of the H atom are the same as 
<2Ial [21bt with W = 0. 

8 Due to ignoring the asymptotic central symmetry of Hamiltonian at r = 0, spatially 
generalized power series if = Y2\ u r ^ u \^{ r \ r 3i • ■ ■ i r Jv)^A/j(^i f) was used by former authors. 
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Therefore the many-electron wave function can be decomposed in the vicinity 
of Coulomb singularities as 



4- = r l [u/(r|r 3j . . . , r N )Y im (i}, if) + r 3 v(r\r 3 , . . . , r N )] , (23) 



where 



ui = u e (0) (l + ar + br 2 + cr 3 + . . .) , 

OO A 

v = E E r X - e - 1 v XlM (r\^...,r N )Y Xll ^,<p). 

X=l+1 fj,=—X 

Term O (r e ) of decomposition ([23)1 is an eigenfunction of the spherically aver- 
aged Hamiltonian, term O (r^ +3 ) reflects anisotropy of molecular or crystalline 
potential, so group theoretical considerations apply only to v = u(r|r3, . . . ,rjv). 

3.4 Local behaviour at Coulomb singularities of spheri- 
cally averaged Hamiltonian 

Let us investigate the local behaviour of first term of decomposition f|23[) in more 
detail. In view of Eq. (|17l) the limiting form of spherically averaged Hamiltonian 
is equivalent to that of a Coulomb potential embedded in a uniform background 

hence Eq. (|19l) can be rewritten as 

ru" + {21 + 2) v! - (2a + ^r) u = 0, 

where a = Mqyq^ and /3 2 = 2M (Wo — E). Seeking the solution in the form u = 
e~P r w(r) and then by making change of variable x = 2/3r we obtain following 
confluent hypergeometric equation 

xw" + {2t + 2 - x) w' - (£ + 1 + a/0) w = 

which is of Kummer type [45] 

xw" + (b — x) w' — aw = 

with a = £+ l + a I f3 and b = 2^+2. Regular solution of this equation is the 
following Kummer function [f| 

where (a) k denotes the Pochhammer symbol defined by (o) fc = T (o+ k) /T (a). 
Therefore the physical solution of Eq. (|19D in the neighbourhood of r = is 

Ui (r) = e-^xFt (£ + 1 + a/f3; 21 + 2; 2/3r) . (25) 



9 The solution reduces to associated Laguerre polynomial in case of Wo = 0. We may 
formally consider —a//3 as a non-integer principal quantum number. 
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Since Eq. (|19[) is non-singular, its analytic solution can be expanded in 
terms of a power series. By substituting power series ut(r) = X)fc°=o a k rk into 
the equation and collecting terms with equal powers of r we arrive at algebraic 
equation 

[(2£ + 2) ai - 2aa ] r° 

OO 

+ H k + : )( 2£ + 2 + fc ) flfe +i " 2aa * " Ak-l] rk = 

leading to recurrence relations 

ai = £^7Y a o> ( 26a ) 
2aa k + /3 2 a fe -i 

a/c+i = i w r, re = 1,2,.... (26b 

+ (2£ + 2 + k)(k + 1)' ' ' V y 

The first one is equivalent to Lowdin's |14j first order cusp condition (TT]) and 
the second one provides a simple way of generating higher order cusp conditions 
for eigenfunctions of the spherically averaged Hamiltonian (|24|) . Use of these 
recurrence relations is more efficient numerically than evaluating the confluent 
hypergeometric function in Eq. (|25l) . 

In view of Eqs. (jlOa[) , (|18p and (l2~5j) the cigenf unction of the spherically 
averaged Hamiltonian is of the form 

il> = r'u e (0)e- f3r 1 F 1 [i + 1 + a/p; 21 + 2; 2/3r) Y lm {d, <p) (27) 

in the vicinity of Coulomb singularities (for comparison with an accurate Hartree- 
Fock-Roothaan atomic wave function see Figures [5] and [3]) . This analytic func- 
tion fully characterizes cusps of spherically symmetric or spherically averaged 
systems since arbitrarily high order cusp relations can be derived from it. Higher 
order cusp conditions obtained from Kummer type wave function ([27]) agree 
with that of Refs. [27 and [30 (after typographic correction Z a i-> Z^)- For 
anisotropic systems, only first and second order cusp conditions (I21al I21b|) are 
exact which can also be obtained from recurrence relations (|26al I26b[) . This 
bound-state local solution becomes non-physical at larger distances satisfying 
Q1Q2/T + Wo > since non-coalescent particles within the sphere of radius r are 
not included. Function (|27p is not square-integrable but knowledge of its higher 
order derivatives at r = is useful in numerical calculation of first term of Eq. 

(E3D- 



4 Many-electron cusp and boundary conditions 



Using expansions (|13[) and (|23p of many-electron wave function about singu- 
lar points of Hamiltonian ([7]) we are able to recover explicit forms of natural 
boundary conditions for the many-electron Schrodinger equation. Due to local 
isotropy ([Sal I9bjl of Hamiltonian at singular points we can follow the same pro- 
cedure as we did in section [2] for the H atom. One of our resulted boundary 
conditions may be considered as an exact cusp condition without spherical aver- 
age which contradicts Kato's result. Using our Kummer function solution (|27l) 
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Figure 2: Wave function of Kummer type (|27l) compared with Hartree-Fock- 
Roothaan (HFR) wave function computed by Koga, Kanayama, Watanabe and 
Thakkar (46] for ground state of the He atom. Radius r — Z/Wo of bound-state 
region of spherically averaged Hamiltonian (|24[) is marked on the r axis. 





Kummer function 




HFR by Koga et al. 




- ,-S 5691 











0.2 



0.4 0.6 

r[a.u.] 



Figure 3: Radial density obtained from (|27p compared with that of computed 
by Koga et al. [32] for ground state of the He atom. Our near-nucleus approx- 
imation is accurate up to surprisingly large distances. Effective Bohr radius r$ 
of the orbital is marked on the r axis. 
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for the spherically averaged Hamiltonian we can define arbitrarily high order 
cusp conditions in the framework of central field approximation. 



4.1 Boundary conditions 

In view of expansion fj 13[) of many-electron wave function at large distances, its 
logarithmic derivative with respect to r is 



y/ Or r— >oo 



M'(Q + 1) 1 



y/-2M'E 

which defines our first boundary condition: 

1 a* 



lim — — = -\f-2M'E. (28) 
i — >oo w (jr 

The directional logarithmic derivative of the many-electron wave function in 
arbitrary direction e can be expressed in terms of spherical polar coordinates as 

e ■ V* e • e r <9* e • e# 5* e ■ e v 9* 

* * dr $r 3t) *r sin dip ' 

where {e r , e#, e v } is basis of the coordinate system. In radial directions, defined 
by e • e r = ±1 and e • e$ = e • e v = 0, above expression reduces to 

e • V* _ e • V r * _ ±1 
* * * dr ' 

Since wave function = r g u(r\T^, . . . , r^v) = r e u(r) given by Eq. (|2"2"|) has a root 
of multiplicity I at Coulomb singularities the l'Hospital rule should be applied 
I times in order to obtain a definite limit 

1 9* d e r +1 ^> I + 1 du(r) , , 
► — > — + O (r) 

dr r- >-0 r->o u(r) dr 



where differentiation rules (|36al I36b[) are used o. Since expansion (123)) exhibits 
hydrogenic angular dependence r ut{r)Yi m { , d : ip) in the vicinity of Coulomb sin- 
gularities the radial logarithmic derivative of "J has the following one-sided lim- 
its: 

1 d^ d e+1 ^ 
lim -^-= lim = ±(£ + l)a, (29) 



where a is defined by Eq. (|21a[) and 

r^+0 = (r-tO, 0,tp), 

r — > —0 = (r — >• 0, 7r — i?, 7r + 

Equation (|29p is our second boundary condition. 

Equations (|2"5|) and ([2l?|) represent homogeneous Robin boundary conditions 
for the many-electron wave function and have the same form as Eqs. ((SJ and 
(|6]) obtained for the H atom which is a consequence of asymptotic isotropy (|9al 
I9b[) of Hamiltonian at singular points. 



10 In general, radial partial derivative of an anisotropic function is anisotropic which ex- 
plains appearance of spherical average in Kato's cusp condition J2{. Sec discussion after the 
statement of theorem I la in Ref. I18| . 
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4.2 Cusp conditions 



Boundary condition (|29l) is the exact form of first order cusp condition which 
does not contain spherical average. By substituting specific values of M, q\ and 
<72 into Eq. (l29|) we obtain cusp conditions for the electron-nucleus coalescence 

^oV^^' ^ = 0,1,2,..., (30) 

where M v — Z v m p /m e + (A u — Z v )m n /m e denotes mass of the ^-th nucleus 
(in a. u.). Similarly we obtain following cusp conditions for electron-electron 
coalescences 

^W = 4 ^0,2,4,..., (31a) 

^o^ff = ± 2' < = W,..., ( 31b ) 
where and ft denote relative singlet and triplet spin states of coalescent 



electrons, respectively. First order many-body cusp conditions (1301 131al I31bp 
are obviously CP-invariant similarly to Eq. ([5]) obtained for the H atom. 

Since second directional derivative (e • V) 2 is of definite sign the second order 
cusp condition has the form 

lim —-7—5- = lim ^— = {£ + 1 ){£ + 2 )b, (32) 

where differentiation rules Q36al I36cl) arc used, b is defined by Eq. (|21b[) . 

Coefficient a in Eqs. ([2U1 13"2"j) depends only on charges and masses of the 
coalescent particles. In view of (TIE)) and (|21b>[) . coefficient b = b (r^, . . . , rjv) in 
Eq. (|32p depends only on charges and distances of non-coalescent particles hence 
it has the same value for the spherically averaged Hamiltonian. Arbitrarily high 
order cusp conditions can be obtained from Eq. (|27p by means of recurrence 
relation (|26a[ I2f3b|) in the framework of central field approximation. Third and 
higher order cusp conditions for exact wave functions require spherical average 
similarly to Kato's cusp condition @ since higher order coefficients of expansion 
(|2"3"]) depend both on distances and directions of the non-coalescent particles. 

In view of Eq. f|23[) the electron pair density exhibits short-range hydrogenic 
angular dependence in the vicinity of coalescence points hence satisfies exact 
first and second order cusp conditions similar to (|31a[ I31b|) and (|32p . 



5 Consequences 

5.1 Physical consequences 

Since jump of logarithmic derivative of wave function is 2Z at the nuclei both in 
the one-electron and in the many-electron cases, its discontinuity is caused by 
common terms of the two Hamiltonians, namely, by singular Coulomb potential 
and by singular kinetic energy of opposite sign as if no other particles were 
present except the coalescent ones. An electron-electron potential taking part 
in cancellation of the nuclear Coulomb singularity would require a more singular 
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wave function than the square-integrable functions hence cusp condition shows 
an evidence for the normalization condition. 

Potential energy term Wo defined by Eq. (|16p is related to the NMR chemical 
shift since it is proportional to ratio Hi (0) /H of induced diamagnetic shielding 
field to the applied external magnetic field [57J EH] . 

Expansion (|23l) characterizes the short-range behaviour of many-electron 
wave function both in the vicinity of electron-nucleus and electron-electron 
coalescences hence it can be used to describe short-range electron correlation 
which depends on the chemical environment in second order through parameter 
b = b(r 3 , . . .,r N ). 

First order cusp condition (|21a[ |2U)) is responsible for boundedness of many- 
electron wave function at Coulomb singularities by exactly cancelling singular 
kinetic and potential energy terms. Second order cusp conditions (|21bl 1321) 
enforce correct value E of local energy Ei oc = ^~ l H^ at Coulomb singularities 
by imposing constraints on curvature of many-electron wave function. These 
conditions are exact for arbitrary values of t and do not need spherical average 
in contrast to Kato's cusp condition ©. 

One can see from Eq. (|27p that all radial partial derivatives of eigenfunctions 
of spherically averaged Hamiltonian at coalescence points depend on the same 
constant quantities: £, qi, 52, M, E and Wq. 

Third and higher derivatives of exact eigenfunctions depend on anisotropic 
multipole potential terms of expansion (Tl5j) as well. In view of Eq. (|2"3"j) . sym- 
metry considerations apply only to terms O (r^ +3 ). 

The real component of Bloch functions exhibits cusps at the nuclei similarly 
to molecular wave functions since the cosine function does not vanish at r = 0. 

5.2 Numerical consequences 

If both boundary conditions are known the numerical solution of Schrodinger 
equation over a grid yields an algebraic eigenvalue problem which can be solved 
by means of Jacobi - Goldstine - Murray - von Neumann diagonalization algo- 
rithm [49] , Before Kato's cusp condition only long-range behaviour (f2"5|) of wave 
function was known hence a numerical trick called shooting method was used 
to substitute the missing boundary condition at Coulomb singularities. This 
method is best suited to equidistant grids, where inward and outward numer- 
ical integrations are equally accurate. Since most of energy is concentrated at 
near-nucleus regions the practical grids used in quantum chemistry are substan- 
tially finer in the vicinity of nuclei than in the interstitial and exterior regions. 
The grid is coarsest at large distances representing 00 in the numerical cal- 
culation hence an outward integrated solution starting from a guessed initial 
condition is fitted at a midpoint to an inward integrated solution based on an 
inaccurate initial condition. Kato's cusp condition ([2]) is suitable only to the 
central-field approximation and cannot be used as an exact boundary condi- 
tion. Our boundary conditions (|28l) and (f2U|) are exact regardless of molecular 
or crystalline symmetry since the wave function exhibits hydrogenic angular 
dependence in both limiting cases. A shooting method based on our cusp con- 
dition (|29p would be more accurate than the traditional one due to the finer 
grid at nuclei. Since both boundary conditions are known one can transform 
Schrodinger equation to an algebraic eigenvalue problem instead of performing 
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the time-consuming shooting loop. Approximate wave function (|27p can be used 
to find optimal near-nucleus step size for the grid. 

Basis sets satisfying cusp conditions improve convergence of Hartrce-Fock- 
Roothaan variational calculations. It is well known that Slater-type basis sets 
satisfy first order cusp condition if exponential decay parameter of one basis 
function is appropriately fixed 

L 

R = c r e e-^ + Y, cxr e+x e-^ r 

A=l 

letting only its weight factor cq to be varied. Similar basis set was used by 
Roothaan and Kelly [50 however their summation inexplicably starts from A = 3 
resulting in a slow convergence for I > hence use of cusp condition was limited 
to the s orbitals in their atomic calculations F*l It is easy to construct a Slater- 
type basis set satisfying both first and second order cusp conditions by equating 
three leading expansion coefficients of linear combination of two Slater functions 
with equal exponential factors to that of expansion (|20p yielding 

R = cor' 

where a and b are defined by Eqs. (|21al [21b]) . In order to preserve asymptotic 
behaviour of wave function as r — > oo a Slater function of the form (|12[) should 
be added to the basis set. 

It is widely believed that Gaussian basis sets are not suitable to describe 
nuclear cusps since they have zero gradients at the nuclei which is true only 
for individual Gauss functions but not for their linear combinations By equat- 
ing three leading expansion coefficients of linear combination of three Gauss 
functions with equal exponential factors to that of expansion (|20j) we obtain 
Gaussian basis set 

4>(x,y,z) = c r e [l + ar + (b + g )r 2 ]e-3° r2 Y em ($,<p) 

L 

+ z k e~^ r2 (34) 

t+j+k>£+3 

satisfying both first and second order cusp conditions, where a and b are defined 
by ([2Tall2Tb| and r 2 = x 2 +y 2 +z 2 . There is no finite linear combination of Gauss 
functions which exhibits asymptotic behaviour (I12[) of wave function as r — > oo. 
Use of basis functions satisfying both first and second order electron-electron 
cusp conditions provides the simplest way to include short-range correlation 
effects. Above basis set is more efficient numerically and requires less modifi- 
cation of existing Gaussian computer codes than implementing a Jastrow-type 
correlation |51] , 

Asymptotic hydrogenic angular dependence of three leading terms of ex- 
pansion (|23[) of many-electron wave function about Coulomb singularities ex- 
plains the success of central-field approximation used in atomic calculations and 

11 We suppose that their summation starting from A = 3 arises from some unpublished 
expansion similar to our Eq. (123 D . where our U( was simply replaced by e~ Zr /(?+ 1 ). 



I + I h — | 



(33) 
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muffin-tin approximation of solid state physics. The Kummer-type confluent hy- 
pergeometric function (|27p intended to characterize the short-range behaviour of 
eigenfunctions of the spherically averaged Hamiltonian is surprisingly accurate 
even at relatively large distances, e.g. one can see from Figure [3] that relative 
error of this approximation to radial density Anr 2 \ip\ 2 is 5.8% at the effective 
Bohr radius and is 0.4% at its half. Therefore the spherically averaged part of 
Hamiltonian is responsible for most of the effects and anisotropic terms can be 
considered as perturbations. 

Appendix: Some differentiation rules 

Leibniz's theorem for differentiation of products states that 
[f(x)g(x)]^ =Y,( n \f^- k \x)g^{x). 



k=0 

For f(x) — x e one obtains 



k (l-n+kV. 

Specific higher order derivatives of the above type used in this paper are 

[x l g{x)] {e) -^J\ g (Q) + 0{x), (35a) 
[x e g(x)f +1) -^U + l)!g'(0) + O(x), (35b) 



l , M(^+2) 



(XT ' — >(* + 2)! ff "(0) + O(a:). (35c) 



For higher order radial partial derivatives of r l u{r) we obtain similarly 

()r '" ,rl ^!u(0) + O(r), (36a) 



fin 

+ O (r) , (36b) 



d e+1 r e u(r) fe | du 



., , , ^+1). ., 

d t+2 r l u{r) d 2 u 
dr i+2 r ^o [ >' dr 2 



r=0 

+ O (r) . (36c) 



Quotients of above derivatives become isotropic if u(r) can be written as a 
product of radial and angular parts. 
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